home *** CD-ROM | disk | FTP | other *** search
/ Practical Algorithms for Image Analysis / Practical Algorithms for Image Analysis.iso / CH_6.1 / XVOR / XVOR.C < prev    next >
C/C++ Source or Header  |  1999-09-11  |  10KB  |  408 lines
  1. /* 
  2.  * xvor.c
  3.  * 
  4.  * Practical Algorithms for Image Analysis
  5.  * 
  6.  * Copyright (c) 1997, 1998, 1999 MLMSoftwareGroup, LLC
  7.  */
  8.  
  9. /*
  10.  * XVOR
  11.  *
  12.  * compute Voronoi diagram or Delaunay triangulation
  13.  * for a set of points in the plane
  14.  *
  15.  * on unsorted data uniformly distributed in the unit square, routine
  16.  * voronoi uses about 20n+140\(srn bytes of storage;  on sorted data, 
  17.  * routine voronoi uses about 160\(srn
  18.  *
  19.  *
  20.  * input: 
  21.  * each input line should consist of two real numbers, separated by 
  22.  * white space.
  23.  * sorted input (option -s):
  24.  * the input is sorted by y coordinate (the second number in the pair);
  25.  * the first input line should be npoints xmin xmax ymin ymax; this
  26.  * describes the number of points and the range of the points; the
  27.  * information is used to determine internal hash table size; it need 
  28.  * not be exact but performance suffers if it is grossly wrong.
  29.  *
  30.  * VORONOI output:
  31.  * there are four output record types.
  32.  * -->s a b indicates that an input point at coordinates a b  was seen;
  33.  * -->l a b c indicates a line with equation ax + by = c.
  34.  * -->v a b indicates a vertex at a b.
  35.  * -->e l v1 v2 s1 s2 indicates a Voronoi segment which is a subsegment of 
  36.  *    line number l with endpoints numbered v1 and v2; if v1 or v2
  37.  *    is -1, the line extends to infinity; e bisects the line segment
  38.  *    delimited by sites s1 and s2
  39.  *
  40.  * DELAUNAY output:
  41.  * output: each output line is a triple i j k, giving the indices of 
  42.  *      the three points in a Delaunay triangle; points are
  43.  *      numbered starting at 0.
  44.  *
  45.  */
  46.  
  47. #include "xvor.h"
  48.  
  49. /* globals */
  50. extern char *optarg;
  51. extern int optind, opterr;
  52.  
  53.  
  54. int VOR = ON;                   /* default: evaluate Voronoi diagram */
  55. int TRIANGULATE = OFF;
  56. int PRE_SORTED = ON;            /* default: assume y-sorted data in input file */
  57. int DBG = OFF;
  58.  
  59. int WRITE_FLAG = OFF;
  60.  
  61. FILE *fpIn, *fpOut;
  62.  
  63. /*
  64.  * scomp()
  65.  *   DESCRIPTION:
  66.  *     sort sites on y, then x, coord
  67.  *     used with qsort!!
  68.  *   ARGUMENTS:
  69.  *     S1: first sort poiint (Point *)
  70.  *     S2: second sort poiint (Point *)
  71.  *   RETURN VALUE:
  72.  *     1, 0, -1 depending on x,y position
  73.  */
  74.  
  75. int
  76. scomp (S1, S2)
  77.      const void *S1, *S2;
  78. {
  79.   const struct Point *s1 = S1;
  80.   const struct Point *s2 = S2;
  81.  
  82.   if (s1->y < s2->y)
  83.     return (-1);
  84.   if (s1->y > s2->y)
  85.     return (1);
  86.   if (s1->x < s2->x)
  87.     return (-1);
  88.   if (s1->x > s2->x)
  89.     return (1);
  90.  
  91.   return (0);
  92. }
  93.  
  94. /*
  95.  * nextone()
  96.  *   DESCRIPTION:
  97.  *     return a single in-storage site
  98.  *   ARGUMENTS:
  99.  *     none
  100.  *   RETURN VALUE:
  101.  *     next site (struct Site *)
  102.  */
  103.  
  104. struct Site *
  105. nextone ()
  106. {
  107.   struct Site *s;
  108.  
  109.   if (siteidx < nsites) {
  110.     s = &sites[siteidx];
  111.     siteidx += 1;
  112.     return (s);
  113.   }
  114.   else
  115.     return ((struct Site *) NULL);
  116. }
  117.  
  118.  
  119. /*
  120.  * freadsites()
  121.  *   DESCRIPTION:
  122.  *     read all sites, sort, and compute xmin, xmax, ymin, ymax
  123.  *   ARGUMENTS:
  124.  *     file: pointer to open FILE
  125.  *   RETURN VALUE:
  126.  *     none
  127.  */
  128.  
  129. void
  130. freadsites (file)
  131.      FILE *file;
  132. {
  133.   int i;
  134.  
  135. /*
  136.  * skip top line in .vin data file
  137.  */
  138.   fscanf (file, "%d %f %f %f %f", &nsites, &xmin, &xmax, &ymin, &ymax);
  139.  
  140. /*
  141.  * read data
  142.  */
  143.   nsites = 0;
  144.   sites = (struct Site *) myalloc (4000 * sizeof (struct Site));
  145.   while (fscanf (file, "%f %f", &sites[nsites].coord.x, &sites[nsites].coord.y) != EOF) {
  146.     sites[nsites].sitenbr = nsites;
  147.     sites[nsites].refcnt = 0;
  148.     nsites += 1;
  149.     if (nsites % 4000 == 0)
  150.       sites = (struct Site *) realloc (sites, (nsites + 4000) * sizeof (struct Site));
  151.   }
  152.  
  153. /*
  154.  * sort data and find extrema
  155.  */
  156.   qsort (sites, nsites, sizeof (struct Site), scomp);
  157.  
  158.   xmin = sites[0].coord.x;
  159.   xmax = sites[0].coord.x;
  160.   for (i = 1; i < nsites; i += 1) {
  161.     if (sites[i].coord.x < xmin)
  162.       xmin = sites[i].coord.x;
  163.     if (sites[i].coord.x > xmax)
  164.       xmax = sites[i].coord.x;
  165.   }
  166.   ymin = sites[0].coord.y;
  167.   ymax = sites[nsites - 1].coord.y;
  168.  
  169. #ifdef SHOW_DATA
  170.   printf ("\n...data after sorting:\n");
  171.   printf ("...parameters: nsites: %4d", nsites);
  172.   printf ("\n   extrema: %f %f %f %f\n", xmin, xmax, ymin, ymax);
  173.  
  174.   for (i = 0; i < nsites; i += 1) {
  175.     printf ("  x[%d] = %f   y[%d] = %f\n",
  176.             i, (sites + i)->coord.x, i, (sites + i)->coord.y);
  177.   }
  178. #endif
  179. }
  180.  
  181.  
  182. /*
  183.  * readone()
  184.  *   DESCRIPTION:
  185.  *     read one site
  186.  *   ARGUMENTS:
  187.  *     none
  188.  *   RETURN VALUE:
  189.  *     site (struct Site *)
  190.  */
  191.  
  192. struct Site *
  193. readone ()
  194. {
  195.   struct Site *s;
  196.  
  197.   s = (struct Site *) getfree (&sfl);
  198.  
  199.   s->refcnt = 0;
  200.   s->sitenbr = siteidx;
  201.   siteidx += 1;
  202.  
  203.   if (fscanf (fpIn, "%f %f", &(s->coord.x), &(s->coord.y)) == EOF)
  204.     return ((struct Site *) NULL);
  205.  
  206. #ifdef SHOW_DATA
  207.   printf ("\n...READONE: fetching site: x[%d] = %f  y[%d] = %f\n",
  208.           siteidx, s->coord.x, siteidx, s->coord.y);
  209. #endif
  210.  
  211.   return (s);
  212. }
  213.  
  214.  
  215. void
  216. exit_messg (str, code)
  217.      char *str;
  218.      int code;
  219. {
  220.   printf ("\n%s", str);
  221.   exit (code);
  222. }
  223.  
  224.  
  225.  
  226. /*
  227.  * usage of routine
  228.  */
  229. void
  230. usage (char *progname)
  231. {
  232.   progname = last_bs (progname);
  233.   printf ("USAGE: %s infile outimg [-w file] [-s] [-t] [-g] [-L]\n", progname);
  234.   printf ("\n%s constructs Voronoi diagram or Delaunay triangulation\n", progname);
  235.   printf ("for a set of points in the plane using Fortune's algorithm.\n");
  236.   printf ("[Fortune 1984]\n\n");
  237.   printf ("ARGUMENTS:\n");
  238.   printf ("        infile: input filename (.vin) (ASCII)\n");
  239.   printf ("                produced by spp or xah; assumes y-sorted\n");
  240.   printf ("                site coordinates unless option -s invoked.\n");
  241.   printf ("        outimg: output image filename (TIF)\n\n");
  242.   printf ("OPTIONS:\n");
  243.   printf ("   -w file : write output to file ( .vdt) as well as stdout\n");
  244.   printf ("         -s: y-sort data in input file\n");
  245.   printf ("         -t: perform Delaunay triangulation\n");
  246.   printf ("         -g: debug mode\n");
  247.   printf ("         -L: print Software License for this module\n");
  248.   exit (1);
  249. }
  250.  
  251.  
  252. void
  253. main (int argc, char *argv[])
  254. {
  255.   struct Site *(*next) ();      /* next is declared to be a ptr */
  256.   /* to a function which returns a 
  257.    * /* ptr to struct Site */
  258.  
  259.   int ic, is;
  260.   int i_arg;
  261.   char *buf;
  262.  
  263.   Image *imgOut;
  264.  
  265.   static char wbuf[20];
  266.   static char *inp_suffix =
  267.   {".xxx"};                     /* default inp file suffix */
  268.   static char *vin_type =
  269.   {".vin"};                     /* suffix for .vin inp file */
  270.   static char *vsuffix =
  271.   {".vdt"};                     /* suffix for V. outp file name */
  272.   static char *dsuffix =
  273.   {".ddt"};                     /* suffix for D. outp file name */
  274.  
  275.   int reset = 1, no_reset = 0;
  276.  
  277. /* 
  278.  * cmd line options
  279.  */
  280.   static char *optstring = "w:stgL";
  281.  
  282.   buf = argv[1];
  283.  
  284. /*
  285.  * parse command line
  286.  */
  287.   optind = 3;
  288.   opterr = ON;                  /* give error messages */
  289.  
  290.   if (argc < 3)
  291.     usage (argv[0]);
  292.  
  293.   while ((i_arg = getopt (argc, argv, optstring)) != EOF) {
  294.     switch (i_arg) {
  295.     case 'w':
  296.       printf ("...option %c: ", i_arg);
  297.       strcpy (wbuf, optarg);
  298.       if ((fpOut = fopen (wbuf, "w")) == NULL) {
  299.         printf ("\n...cannot open %s for writing\n",
  300.                 wbuf);
  301.         exit (1);
  302.       }
  303.       printf ("write data to file %s\n", wbuf);
  304.       WRITE_FLAG = ON;
  305.       break;
  306.     case 's':
  307.       printf ("...option %c: ", i_arg);
  308.       printf ("perform sort on y-coord. of input\n");
  309.       PRE_SORTED = OFF;
  310.       break;
  311.     case 't':
  312.       printf ("...option %c: perform Delaunay triangulation\n", i_arg);
  313.       VOR = OFF;
  314.       TRIANGULATE = ON;
  315.  
  316.       /* construct new output file name */
  317.       ic = 0;
  318.       while ((*(wbuf + ic) = *(buf + ic)) != '.')
  319.         ic++;
  320.       for (is = 0; is < 4; is++)
  321.         *(wbuf + (ic) + is) = *(dsuffix + is);
  322.       break;
  323.     case 'g':
  324.       printf ("...option %c: debug mode\n", i_arg);
  325.       DBG = ON;
  326.       break;
  327.     case 'L':
  328.       print_sos_lic ();
  329.       exit (0);
  330.     default:
  331.       printf ("...unknown condition encountered\n");
  332.       exit (1);
  333.       break;
  334.     }
  335.   }
  336.  
  337.   printf ("read data from file %s\n", buf);
  338.   /* get size of data set from file */
  339.   if ((fpIn = fopen (buf, "r")) == NULL) {
  340.     printf ("\n...cannot open input file %s\n", buf);
  341.     exit (1);
  342.   }
  343.  
  344.   /* construct output file name */
  345.   ic = 0;
  346.   while ((*(wbuf + ic) = *(buf + ic)) != '.')
  347.     ic++;
  348.   for (is = 0; is < 4; is++)
  349.     *(wbuf + (ic) + is) = *(vsuffix + is);
  350.   /* strip input file suffix */
  351.   for (is = 0; is < 4; is++)
  352.     *(inp_suffix + is) = *(buf + ic + is);
  353.  
  354. #ifdef PORTING
  355.   printf ("\n...give size of all structures defined in header file:\n");
  356.   printf ("   struct Freenode: %d\n", sizeof (struct Freenode));
  357.   printf ("   struct Freelist: %d\n", sizeof (struct Freelist));
  358.   printf ("      struct Point: %d\n", sizeof (struct Point));
  359.   printf ("       struct Site: %d\n", sizeof (struct Site));
  360.   printf ("       struct Edge: %d\n", sizeof (struct Edge));
  361.   printf ("   struct Halfedge: %d\n", sizeof (struct Halfedge));
  362. #endif
  363.  
  364. /* 
  365.  * get data from input file
  366.  */
  367.   if (strcmp (inp_suffix, vin_type) == 0) {
  368.     freeinit (&sfl, sizeof (struct Site));
  369.     if (PRE_SORTED == ON) {
  370.       fscanf (fpIn, "%d %f %f %f %f",
  371.               &nsites, &xmin, &xmax, &ymin, &ymax);
  372.  
  373.       next = readone;
  374.     }
  375.     else {
  376.       freadsites (fpIn);
  377.       next = nextone;
  378.     }
  379.  
  380.     printf ("\n...initialize geometry settings...\n");
  381.     siteidx = 0;
  382.     geominit ();
  383.  
  384.     /* allocate default image buffer for image output */
  385.     imgOut = ImageAlloc ((long) (ymin + ymax), (long) (xmin + xmax), BPS);
  386.  
  387.     printf ("\n...initialize plotting...\n");
  388.     plotinit (imgOut);
  389.  
  390.     if (VOR == ON)
  391.       printf ("\n...start voronoi construction...\n");
  392.     else if (TRIANGULATE == ON)
  393.       printf ("\n...start Delaunay triangulation...\n");
  394.     voronoi (next, imgOut, WHITE);
  395.  
  396.   }
  397.   else
  398.     exit_messg ("\n...unknown inp file suffix", 1);
  399.  
  400.  
  401.   if (WRITE_FLAG == ON)
  402.     fclose (fpOut);
  403.   printf ("\n...writing graphics output file %s\n", argv[2]);
  404.   ImageOut (argv[2], imgOut);
  405.  
  406.   fclose (fpIn);
  407. }
  408.